NASA-CR-190667 


/ 


OSU-ECE Report NASA 92-02 


/ * - 
/ / ^ 



Phase 2 


Annual Report on 


NONLINEAR STABILITY AND CONTROL STUDY 


OF HIGHLY MANEUVERABLE 
HIGH PERFORMANCE AIRCRAFT 


rsj 

N 

cr 

0 

1 

rg 




U 

C 


(NASA Grant No. NAG-1-1081) 
Date: August 20, 1992 
R.R. Mohler, Principal Investigator 


Oregon State University 

Department of Electrical and Computer Engineering 
Corvallis, Oregon 97331-3211 

(503) 737-3470/3617 


c 

u. 3 

3 rvi 

C' 


* > xi 
< O t-o 
UJ 3 < 

y- x x 

h ’/* O a. 
_j *— 

-J X - 

C' l- 

Z ^ : .J UL 

-J «1 


<r 


c 


1 


c.j C/ 



r~i 3 0 C 

i > t rr - 

c£ > < - v 


t -J >- c. *— 

<t -j 3 r 

t/“' r ~T Li. "D C * 
«.< < 3 X C 

3 * :i .* C 

*) X .X * s' 


Graduate Research Assistants: S. Cho, C. Koo, R. Zakrzewski 
Undergraduate Participants (NSF Support): D^Aaberg, J. Young 
Visiting Researchers: J. Dory*, J. Kurek , A. Yagen 
NASA Fellow: D. Collins 


*ADA-lsrael Support 
^International Exchange Board 


G3/0B 0115^70 


OSU-ECE Report NASA 92-02 


Phase 2 

Annual Report on 


NONLINEAR STABILITY AND CONTROL STUDY 
OF HIGHLY MANEUVERABLE 
HIGH PERFORMANCE AIRCRAFT 


(NASA Grant No. NAG-1-1081) 
Date: August 20, 1992 
R.R. Mohler, Principal Investigator 


Oregon State University 

Department of Electrical and Computer Engineering 
Corvallis, Oregon 97331-3211 

(503) 737-3470/3617 


Graduate Research Assistants: S. Cho, C. Koo, R. Zakrzewski 
Undergraduate Participants (NSF Support): D^Aaberg, J. Young 
Visiting Researchers: J. Dory*, J. Kurek**, A. Yagen 
NASA Fellow: D. Collins 


*ADA-Israel Support 
^International Exchange Board 


TABLE OF CONTENTS 


Page 

1. OVERVIEW 1 

2. NONLINEAR MAC ALGORITHM STATUS 7 

3. TIME-OPTIMAL CONTROL 10 

3.1 Introduction 10 

3.2 Switching-Time-Variation Method H 

3.3 Approximation for Systems Affine in Control 14 

3.4 Computer Implementation 15 

4. SLIDING-MODE CONTROL AND FEEDBACK LINEARIZATION 19 

4.1 Sliding-Mode Control 22 

4.2 Feedback Linearization 22 

4.3 Continuous Sliding-Mode Control with Feedback Linearization 26 

4.4 Simulation 30 

5. OTHER PRELIMINARY STUDIES 37 

5.1 Bilinear-Based Adaptive Control 37 

5.2 Neural-Net-Based Control 38 

6. REFERENCES 41 

APPENDICES 

A. List of Project Publications 

B. On Nonlinear Model Algorithm Controller Design 

C. Analysis of Nonlinear Stability Using Robust Stability Analysis for Linear Systems 


i 


1. OVERVIEW 


This research should lead to the development of new nonlinear methodologies for the adaptive 
control and stability analysis of high angle-of-attack aircraft such as the FI 8 (HARV). The present 
progress report reviews project research over the second year. 

The emphasis has been on nonlinear adaptive control, but associated model development, system 
identification, stability analysis and simulation is performed in some detail as well. Table 1 summarizes 
various models under investigation for different purposes. 

Models and simulations for the longitudinal dynamics have been developed for all types except 7 
in Table 1. A very preliminary analysis has been made on type 6 (neural net models) for adaptive control 
thus far. It has been shown that dynamic accuracy roughly increases with ascending order of model type 
from 1 to 7, except that perhaps 3 (Volterra series) and 6 (neural nets) should be interchanged. However, 
such comparisons depend on how the models are utilized. For example, the neural-net model and 
subsequent control seems to be more accurate than 1 to 5 for flight profiles outside & priori data sets. 
Here, the focus is on adaptive control, generated by model -reference types 1 to 6, of a complex nonlinear 
aircraft motion represented by 7 (nonlinear ordinary differential equations). Preliminary analyses use a 
nonlinear second-order approximation [1] which we found useful for changes in angle of attack (a) by 
about 10° or possibly 20°. A fifth-order nonlinear longitudinal model with the traditional stability 
derivatives generated as functions primarily of a, for a given altitude and mach number, successfully 
mimicked F18 flight trajectories [2], and is being utilized for our nonlinear adaptive-control studies at 
the present time. These models are discussed in the project’s first annual report [3]. Simulations and 
studies reported here refer to the complex (fifth-order) longitudinal model unless otherwise noted. 

Briefly, studies completed indicate that nonlinear adaptive control can outperform linear adaptive 
control for rapid maneuvers with large changes in a. Figures 1 and 2 compare the transient responses 
where the desired a varies from 5° to 60° to 30° and back to 5° all in about 16 sec. Here, the 
horizontal stabilator is the only control used with an assumed first-order linear actuator with a 1/30 sec 
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Table 1. Aircraft Models 


Linear perturbations de- 
sired a, M, h 


Gain scheduled (nonlinear 
function of a) from 1 


Volterra seriesf 

a) at reference states 

b) general case 


Bilinear system 

a) continuous 

b) BARM A 

Polynomial time series 


Neural network 


Nonlinear ordinary differ- 
ential equation model 


Nonlinear partial differ- 
ential equations model 


Purpose 

Remarks/Limitation s 

Local control, check of nonlinear system, 
application of well developed linear con- 
trol methodologies 

Local stability 

Only valid for small maneuvers 
Special case of types 2-5 

Gain-scheduled adaptive control based on 
well developed methodologies 

Simplified description of complex system 

Approximate stability 

May have stability problems with 
small number of reference states 
and/or large fast maneuvers 

Nonlinear adaptive control via cross-corre- 
lation and/or & priori dynamic structure 

Stability approximation 

Simplified dynamic description of complex 
system 

Non-orthogonal series approximation 

Sufficiency of 2 or 3 kernels 

Large computation time for adapta- 
tion 

Nonlinear adaptive control via model 
reference identification (NLMRAC) 

Stability approximation 

Simplified dynamic description 

Large computation time 

Bilinearizing controllers may be 
more practical than linearizing ones 

Polynomial approximation may be 
more accurate but more time con- 
suming than linear or bilinear ap- 
proximation 

Potential application to adaptive control 

Probably less accurate than 4 or 5 
for a given data set but accuracy 
may be more robust outside the 
available data set 

Accurate approximation to fast large ma- 
neuvers for "final" design and simulation 

Stability 

Neglects flexible modes and other 
complications 

Allows the treatment of distributed exter- 
nal stresses and local pressures as well as 
internal stresses and deformations 

Computations are time consuming 

Overall model complexity 

Allows one to treat flutter and dis- 
tributed control strategies, etc. 


is angle of attack, M mach number, and h altitude 
f Wiener series can be used for orthogonal representation 
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Figure 1. Response for Linear MAC 
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Figure 2. Response for Nonlinear MAC 
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time constant. Unfortunately, an additional rate constraint significantly reduces the system performance 
for both the nonlinear and linear adaptive control as shown in Figures 3 to 5 and analyzed in the next 
section. Such lack of controllability can be improved, of course, by introducing thrust vectoring as used 
in [2]. Appropriate thrust-vector control to supplement the traditional pitch-motion stabilators is 
underway for the nonlinear adaptive controllers, and preliminary results are encouraging. 

A preliminary analysis of time-optimal control of a is studied in Section 3. Here, a new algorithm 
is derived from the switching-time variational method [4,5] and then applied successfully to the simplified 
second-order nonlinear model [1]. The method has been adapted to the more complex nonlinear fifth- 
order model and new simulations are in progress. This study should provide a "yard stick" by which to 
evaluate controller performance as well as provide a base for more effective controller designs. For 
example, it might be used in conjunction with a neural-net generated approximately, time-optimal 
feedback controller as introduced in Section 5.2. As a byproduct of this analysis the complicated 
Jacobian of the longitudinal dynamics will be computed as a function of a and other variables. While 
it is used here to compute bang-bang controller switching times, it should have other uses for approximate 
dynamic-system identification beyond the usual time-invariant linearized models at trim states. 

Nonlinear PIF (proportional plus integral plus filter) control as a generalization of the linear PIF 
control is being studied. The linear PIF controller developed by Ostroff at NASA-LaRC has shown 
remarkable success in simulations as demonstrated in [4] and by more recent simulations and flight tests 
at NASA-LaRC on HARV. The first nonlinear PIF controller which we are investigating involves a 
bilinear (linear with parametric control) model reference for basic control philosophy. For example, the 
PIF model reference may be established from a simplified/approximated time-optimal solution. Then a 
bilinearization of the complex nonlinear aircraft about say an approximate time-optimal trajectory can be 
used to generate the PIF control by possibly an optimal control such as with a quadratic performance 
index. This is introduced briefly in Section 5.1. 
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Figure 3. Response of Linear MAC with Stabilator Limit of 40° per Second 


nonlinear MAC response limiting stabilator rate to 40 deg/sec 
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Figure 4. Response of Nonlinear MAC with Stabilator Limit of 40° per Second 
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Figure 5. Response for Nonlinear MAC with Stabilator Limit of 40° per Second and Slower Reference 
Trajectory 

Feedback linearization and variable-structure or sliding-mode control also are being investigated 
in this project. Both generally have major robustness problems for complex nonlinear dynamics with 
major unknown components and control constraints such as HARV due to aerodynamic stability-derivative 
uncertainty and available actuators. However, they might be useful as a foundation for more complex 
nonlinear controller design synthesis. Also, they can be useful to provide benchmarks of performance 
for comparison with time-optimal controls and other nonlinear controls studied here. At the other 
extreme, a comparison with robust linear controllers will also be made. 

A cursory analysis (presented in Section 4) suggests that robust linearizing controllers (as well as 
robust linear controllers) of sufficient complexity may very well stabilize the longitudinal HARV motion 
for even large deviations from trim flight, but there most likely would be considerable reduction in 
maneuverability along with transient overshoot — the two of which exhibit some tradeoff. It should be 


control signal (dashed) and stabilator angle 



10 15 20 25 30 

time (seconds) 




6 


noted here that linear adaptive control really involves nonlinear feedback which can be thought of as a 
bilinear system (linear system with parametric control). In conjunction with linear model -reference 
adaptive control, the nonlinear feedback essentially becomes a feedback-linearization process, and again 
may suffer the corresponding drawbacks in practice. 

Also, it should be noted that ideally, at least, the sliding-mode control (with probable chatter) can 
conveniently be implemented to accommodate state constraints (as an alternative to excessive overshoot) 
such as can arise from the time-optimal case. But, again, its feasibility with the physical actuator’s rate 
and accelerating constraints needs to be studied further. 

Finally, linear stability arguments are developed in Appendix C which tend to at least an 
approximation of the admissible range of model parameters as applied to the nonlinear second-order 
approximation [1]. 


2. NONLINEAR MAC ALGORITHM STATUS 

Model algorithmic control (MAC), described in [3], starts with 

«ref< k+1 ) = a mod( k + 1 ) + («( k ) ' a mod( k >) 

where 


“mod( k + 1 ) = P a T< « k ) 


( 2 ) 


m = 



(3) 


As the control at the moment k must be already computed at moment k the values of a(k) and q(k) are 
not available for its computation so their estimates must be used instead. The correction term is taken 
to be the prediction error from the moment k-1 and the equation becomes 

««f(k + l) = “mod( k + D + H k -D - «mod( k -D) (4) 
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with 


<W k + 1 ) = 

0(k) = [d,a 2 ,a 3 ,q,qa,qa 2 ,qd 3 ,u,ua,ua 2 ,uor 3 ,l] (k) 
a(k) = pj<^(k-l) + (a(k-l) - a mod (k-l)) 

m = p q T ^>( k -i) + (q( k -D - q m od(k-D) 

The controller is assumed to know the values of angle of attack and of pitch rate at the moment k-1. 
Then it estimates their current values a(k) and q(k) taking into consideration previous prediction errors, 
and based on them it calculates the control required to achieve a ref at the moment k+ 1. The value of 
control is found as: 

... d r - p la or - p 2a & 2 " P 3a « 3 " P4«A - P5a^« ~ Pfe*^ 2 " P7a^“ 3 ' Pl2a ^ 

u(k) ; 72 — 

P8a + P9a“ PlOa« + Plla« 


where 

«r = «ref< k+1 ) “ " “mod^" 1 )) (6) 

and a = &(k), q = q(k) as described above. 

This algorithm was made to be adaptive, or self-tuning, by incorporating on-line identification of 
the parameters. A recursive least squares (RLS) algorithm was implemented in the form taken from 

p(k) = QOc- 2 )^ (7) 

X(k-l) + 0(k-l) T Q(k-2) tf>(k-l) 


Q( k -D 


1 

WT) 


Q(k-2) 


Q(k-2) <ft(k-l) </>(k-l) T Q(k-2) 
X(k-l) + <^>(k-l) T Q(k-2) </>(k-l) 


( 8 ) 


e(k-l) = y(k) - p T <Kk-l) 


( 9 ) 
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where y may denote a or q and p may stand for p a or p q , respectively. The forgetting factor X was 
introduced to enable the algorithm to change the estimates of parameters with the change of operating 
conditions. To avoid the unlimited growth of covariance matrix Q at the steady state when the input is 
not persistently exciting, the variable forgetting factor policy was implemented: 

X( k) = i - £ w 


where e(k) is the current prediction error, e(k) is the average prediction error from last 10 samples, and 
e is equal to 0.01. As an additional precaution, the trace of the covariance matrix Q was monitored and 
Q was reset to diagonal matrix whenever the threshold value was exceeded. 

To further damp the response, the controller is designed to minimize the one step ahead cost 

function: 


J = (y mo d( k+1 ) - y r( k )f + pto ■ “(k- 1 )) 2 
with y mod , y r as before. Minimization of (1 1) with respect to u(k) yields 

(y - a)b + pu(k-l) 
u(k) = L 

b 2 + p 


(ID 


( 12 ) 


where 

a = Pi a a + p 2a « 2 ♦ p 3a ot 3 + p^q + p 5a qa - Pt^" 2 + + Pl2a 

b = Pg a + P 9a “ + Pl0a“ 2 + Pilot® 3 

Obviously, for p = 0 (12) reduces to (5) while for p = « we have u(k) = u(k-l) = const. 

This controller is used in Figures 2, 4, and 5 with only the linear portion of a,b used in Figures 
1 and 3. The algorithm will be generalized to include thrust vector control and variable-horizon cost. 
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3. TIME-OPTIMAL CONTROL 
3.1 Introduction 

Various control strategies have been developed by the team and to find their merits it seems useful 
to have an idea of what are the best output and state trajectories theoretically possible, given the existing 
constraints on the control variables. For substantially nonlinear systems the problem of synthesis of the 
optimal feedback control law is usually untractable. On the other hand, there exist numerical techniques 
that allow us to calculate "open loop controls” - i.e., the specific control signals necessary to achieve the 
minimum performance index. Aware of the difficulties connected with the controller synthesis problem 
we do not seek its exact solution; at this time, we merely want to find the limit for the performance of 
a controller assuming perfect knowledge of plant dynamics and absence of any unforeseen disturbances. 

This report is concerned with the problem of time-optimal control in which we are interested in 
transferring the system’s state from an initial value to some prescribed terminal set in minimal time. In 
the aircraft problem this might mean changing the flight’s pitch angle, path angle, or angle of attack from 
an initial equilibrium value to some other terminal value, preferably also with all other states moving to 
the equilibrium. The control value (stabilator or elevator angle) is naturally bounded from below and 
from above. For some systems it turns out that in case of such simple cube-type constraints on control 
variables, the time-optimal control is of bang-bang type. However, for quite a large class of systems that 
are affine in control, we may approximate any measurable control signal with a bang-bang signal with 
arbitrary accuracy in the sense that corresponding state trajectories are arbitrary close to each other in 
L 1 metric. Hence, also time-optimal control, if it exists, may be approximated by a bang-bang control, 
even if it contains singular arcs. Therefore, the approach presented here is to find the bang-bang control 
that will minimize the transition time. The computational algorithm used here is the switching-time- 
variation method developed in [4,5]. Since the algorithm gives as an output a control signal with finite 
number of switchings, it is tacitly assumed that with large enough finite number of switchings, we are 
able to achieve good enough approximation of optimal control. This, unfortunately, does not follow from 
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the theory I am aware of, since the above mentioned approximation result holds only for bang-bang 
signals with possibly infinite or even uncountable number of switchings. This delicate question is left 
aside for the time being to be clarified later. Another point worth indicating here is that resulting control, 
in an attempt to approximate a continuous "singular" control, may have inter-switching times very small, 
thus precluding any practicality of the approximation. This, however, is of no concern to us since, as 
mentioned before, we are interested only in finding the best possible output, or state, trajectories - not 

the actual control signals corresponding to them at this time. 

A computer program has been developed for numerical solution of the problem. The program, due 
to its modular construction, easily allows various plant models to be plugged into it. The switching time 
variation method is used in it for fixed terminal time with the quality function being the weighted distance 
of the target set. Then the smallest such time is found that allows it to hit the target exactly, and finally 
the optimal number of switchings is iteratively found that gives minimal transition time. 

In what follows the switching-time-variation method is briefly characterized in Section 3.2. Section 
3.3 discusses briefly the approximation theorem for bang-bang controls in systems affine m control. The 
computer implementation of the algorithm is discussed in Section 3.4. Section 3.5 contains the test 
results of the program for a second-order model of longitudinal dynamics of an aircraft. The concluding 
remarks discuss the possibilities of application of the computer package to solutions of more complex and 
problems more close to reality. 

3.2 Switching-Time-Variation Method 

The switching-time-variation method used here was taken from [4], and the original thesis [5] was 
also consulted for the details. The method is designed for the computation of optimal control in the class 
of bang-bang control signals with finite number of switchings. The quality criterion is assumed to be 

J = j^ f (f 0 (x) + g 0 (x) u(t))dt (13) 
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for the system of the form 


— = f(x) + g(x)u(t) 


(14) 


where x e R n , ufR l ,tc [t^tf]. To ensure the existence and uniqueness of solutions of (14), f and g 
are assumed to be continuously differentiable with respect to x. The control values are constrained by 


Of course, any control constraints of the cube-like type u min < u(t) <S Umax may be transformed to form 
(15) for system affine in control. The control objective is to minimize the quality criterion (13) for given 
initial state x 0 with possible penalty term connected with final state already included in f 0 and g 0 by 
standard transformations, assuming that admissible controls are bang-bang with finite number of 
switchings. The version of the algorithm described in [4] was developed for systems with scalar controls 
and the computer program described here is also designed for this special case. However, it is not of 
any particular difficulty to generalize the algorithm to the case of u t R m . If the need arises, the 
computer program may also be modified to accommodate this possibility. Here the scalar version will 
be presented because of its notational simplicity. 

The method is an iterative one - in each step the gradient of the quality criterion with respect to 
switching times is computed. The switching vector is defined as: 


where N is the number of switchings, with constraints. 

to ^ Tj < ... ^ r N < t f 

The control value on the interval [r,,r 1 + 1 ) is then equal (-1) 1 . The augmented system is defined as 


= f(x) + g(x) u(t) 


(18) 
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where F = (f 0 ,f T ), F = (g 0 ,g T ), x£ = (0,xj), and the adjoint system equation is 


dX 

dt 


r |l(t) + ||(t)u(t) 
3x ox 


( 19 ) 


with terminal conditions \(tf) = (dJ/dxjXtf). Then the gradient of the quality criterion with respect to 
the switching vector may be calculated by means of the formula 


with function <f> defined by 

<t>( t) = 2{g(x(t)),x) 


( 20 ) 


( 21 ) 


with gradient calculated the method consists of iterative descent steps 


r ^k-l) 


= * k) + k ‘^ 


( 22 ) 


where kj are such that constraints (17) are satisfied and sufficiently small to ensure that J(k+ 1) < J(k). 
The algorithm is terminated if either the gradient is zero or no feasible (i.e., descent) step may be 

executed. 

On top of the algorithm of finding the optimal switchings with their number given there is an outer 
loop modifying this number. If the optimal control results in a constraints r ; ^ r i + 1 active than the 
switchings i and i+ 1 should be removed. On the other hand, if there are two zeros of <K t) not coinciding 
with any of the switching times than two switchings should be added between these zeros. After the 
modifications of the dimensionality of the switching vector the inner loop of optimization is again 
performed and the process is terminated when no more changes of the number of switches are necessary. 

It is worth noticing that the above algorithm of finding the optimal bang-bang control may be also 
generalized for broader class of systems dx/dt = f(x(t),u(t)). The main difference would be the formula 
for function <£(t0) which would become 
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m = (m > u max) min )) 9 


(23) 


Of course, the technical assumptions ensuring the existence of solutions should be satisfied. 

3.3 Approximation for Systems Affine in Control 

The algorithm described above calculates the optimal control within the class of bang-bang control. 
However, for systems affine in control a result is available stating that we may approximate an arbitrary 
admissible control with a bang-bang control such that corresponding trajectories are arbitrary close. 

The theorem, stated and proven in [6], assumes that we have a system of the form (14) with 
constraints (15). Functions f and g are continuously differentiable, and a Lipshitz type condition 
<f(x) + g(x)u,x > < K(1 + | x || 2 ) preventing finite escape-time is also assumed to be satisfied for all 
x in the region of interest. Then an arbitrary measurable control signal u(t), t e [t^tf] satisfying (15) is 
considered with corresponding state trajectory x(t). Then the theorem states that given any e > 0 is 
always possible to find a bang-bang control u*(t) satisfying |u*(t)| £ 1 , such that the corresponding state 
trajectory x*(t) approximates x(t) uniformly on with accuracy less than e, i.e., |x(t) - x*(t)| <£ e 
for all t e [t^tf]. 

Although the theorem stated above considers a bang-bang control with not necessarily finite or even 
countable number of switchings, it gives some justification to using the switching-time-variation method 
for systems with singular optimal controls. Intuitively for reasonably smooth systems there should be 
some kind of continuity enabl ing in nim approximating the bang-bang control u* with a sequence of bang- 
bang signals u with finite number of switchings. However, we are not aware of any such result, and in 
monograph [7] from 1990 the aforementioned result is cited after [6] as the only one available. It still 
seems feasible to come up with some, maybe more restrictive, assumptions which would justify using 
finite number of switchings. 
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3.4 Computer Implementation 

The algorithm discussed in Section 3.2 was implemented as a quite general software package. It 
finds the time-optimal control for the case when the terminal set is a single point y. The time-optimal 
problem with fixed terminal set is replaced with a sequence of fixed time and free terminal state problems 
with quality index 

J = E Pi (xi(t f ) - y ; ) 2 (24) 

Switching-time-variation method is used to solve this problem, and the desired final time is decreased if 
the resulting quality is zero or is increased in the opposite case. This iteration is repeated until we get 
to the limit time tj below which the quality is always positive, i.e., it is not possible to find a bang-bang 
control transferring the system from Xq to y. 

The optimization method described in Section 3.2 was modified somewhat in details of the gradient 
minimization routine. Instead of performing single step in the direction, a directional search is performed 
with constrained step size. A combination of two-point gradient parabolic approximation and three-point 
non-gradient parabolic approximation is used to find the minimum in the direction. The generation of 
the descent direction is also somewhat different. First, if any of the constraints (15) are active and the 
gradient points outside the feasible region, the gradient is projected on the proper constraining 
hyperplane. The special structure of constraints causes the projection to consist solely of putting the 
appropriate coordinates of the gradient to zero. Then the direction is tangent to the constraining 
hyperplane, and we get an optimization problem of reduced dimensionality. This problem is solved using 
a conjugate gradient method in the version proposed in [6]. The conjugate gradient is restarted not only 
every N iterations, where N is the current dimensionality of the problem, but also whenever the set of 
active constraints changes - i.e., when the algorithm hits or leaves a constraining hyperplane. The 
termination of the procedure occurs when the projected gradient is zero - i.e., no feasible descent step 
is possible, or equivalently when the dimension of the current optimization hyperplane becomes zero. 
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The calculation of the quality criterion and of its gradient involves numerical integration of Eqs. 
(18) and (19). This is done using a fourth-order Runge-Kutta integration method. To integrate the 
adjoint equation (19) the whole state trajectory resulting from integrating (18) must be stored, but for 
calculation only a small number of points from the costate trajectory is needed. 

The program is written in the fashion enabling easy substitutions of different plant models and 
different optimization tasks. To use another model one has simply to provide the routines calculating the 
right-hand sides of Eqs. (18) and (19). The problem is defined in a straightforward fashion by setting 
the values of initial state, terminal state, initial estimate of the final time, etc., in the mam routine. The 
whole program is written in C programming language, and although compiled and run on an IBM PC, 
it may be easily ported to any machine with C language compiler. The only difficulty that may occur 
with more complex systems is the rather severe storage requirements - whole state trajectory has to be 
stored with sufficiently small discretization step in order to calculate the gradient. And, of course, there 
will always be the problem with the speed of calculations for higher dimensional systems. 

The program described above was tested on a model previously used (in our NASA project), i.e., 
a simplified longitudinal aircraft model of second order described in [6]. For the aircraft model the 
problem solved was to increase or decrease the value of angle of attack between trim states. The control 
signal, the elevator angle, was assumed to be between 0 and -20 degrees. A series of maneuvers were 
simulated and reported in [8]. A typical case is shown in Figures 6a and b for maneuvers between 0° 
to 20° and 20° to 0°. It may be observed that the time-optimal control has only one switch for the 
simplified nonlinear model, but the trajectory for a had a substantial overshoot, and consisted of an 
almost linear first portion with high slope before the switch and of slowly decreasing second portion. 

The computer program presented here is suitable for calculating the time-optimal controls for 
arbitrary finite-dimensional systems which are affine in control. The simulation results discussed here 
have mainly testing value showing that the program is in operation. The next step, which is nearing 
completion, uses the program on more complex models such as a fifth-order longitudinal-aircraft model 
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with a linear actuator and affine in control. The resulting time-optimal trajectories for different 
maneuvers could be used as benchmark tests for other controllers or as reference trajectories for time- 
series-based, adaptive, PIF, one-step-ahead (or many-steps-ahead) control. 

For large variations in pitch and angle of attack of the complex model it is found that numerous 
switchings are required and it appears that for extreme cases a singular solution exists for the time- 
optimal control. (Such solutions evolve when extremal trajectories exist with a zero-valued control 
switching function.) To handle such cases more efficiently and more accurately, a new algorithm was 
developed which solves the optimal control problem for a sequence of fixed terminal times with quadratic 
criteria in errors about the desired terminal state. 

Figures 7 and 8 show simulated examples of F-18 maneuvers from trim with pitch 
0(0) = a(0) = 5°, stabilator angle 5(0) = -7°, and constant thrust force of 2956 pounds without thrust 
vectoring and actuator velocity constraints beyond the 30 sec' 1 lag. Figures 7 and 8 show the simulated 
time-optimal maneuvers from this trim state to 6 = 40°, a = 30°, q = 6 = 0 in about 1.4 sec. and to 
6 = 80°, a = 72° in about 1.9 sec. These results are considered to be of a preliminary nature with 
further system complications (thrust vector control, actuator rate limits, etc.) being incorporated into the 
analysis. 

4. SLIDING-MODE CONTROL AND FEEDBACK LINEARIZATION 

Sliding-mode control and feedback linearization can be used in an integrated manner sometimes to 
control effectively nonlinear dynamical plants. Unfortunately, both can suffer considerable performance 
deterioration in the presence of dynamic uncertainty and control constraints for some physical applications 
[3,9]. These aspects need further study with the material presented here only a start. 
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4.1 Sliding-Mode Control (SMC) 

The sliding-mode control, SMC, is a variable high-speed switching feedback control. The main 
advantages of SMC are robustness to plant dynamics and disturbances [10-20]. However, SMC may 
suffer robustness in performance as a result of unknown plant dynamics and control rate or magnitude 
constraints. 

Existence of a sliding mode requires stability of the state trajectory to the sliding surface s(x) = 0 
at least in a neighborhood of {x|s(x) = 0}. The largest such neighborhood is called a region of 
attraction. The existence problem of a sliding mode can be based on the second method of Lyapunov. 
That is, a general Lyapunov function, V(t,x), which is positive definite and has a negative first time 
derivative in the region of attraction. Consequently, the existence and uniqueness of solutions to SMC 
results from the method of Filippov, such that [10-13] 

s(x)s(x) < 0 . 


4.2 Feedback Linearization 

Feedback linearization is an approach to control a nonlinear system by transformation into a fully 
or partly linear one so that linear control techniques can be applied [13,18]. In this respect, consider the 
single input, single output nonlinear system 

* * f < x ' u > • (26) 
y = h(x) 

where x E R n ; u, y are scalars. 

Our objective is to make the output y(t) follow the desired trajectory y d (t) within the bound of the 
state, where y d (t) and its derivatives are assumed to be known and bounded. Here input-output 
linearization, in which relative degree is less than system dimension n, is considered. In this case, it is 
assumed that zero dynamics which is the dynamics of the system subject to the constraint that the output 
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is zero, must be exponentially minimum phase. The main idea is to repeatedly differentiate the output 
y until the output y is related to the input as follows: 

y = Vh(f + gu) = L f h(x) + L g h(x)u @7) 

where Lgh(x) is a Lie-derivative [13,18]. 

If L g h(x) is not equal to zero for all x in the neighborhood x G , we should differentiate until 
L g L[ _1 h(x) ^ 0 for some integer r which is relative degree of r of the system. Then the control law 

u = L (-L f r h(x) + v) (28) 

L g L f r_1 h(x) 


results in the linear differential relation 


y 


(r) 


= v 


(29) 


When r < n, there exists a diffeomorphic coordinate transformation (f,ij) — T(x) for the normal 


form 


1 

I 


I 

N> 

1 

ir-1 

tl 

fr 

i 

-t * 

I 


_a(ft?/) +b(ftjj)u_ 


l = w(ftjj) y = h(x) = ft 


where 


(30) 


(31) 


r = 


h L f h ... L f r_1 h 


V = [ill -.fln-r] 1 
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a(to) = L f r h(x) 


( 32 ) 


b(f,i?) = L g L f r h(x) 


Now the control law is 


u 


1 

b(i\i?) 


[-a(^) + v] 


which results in the linear differential relation 


y 


(r) = 


V 


If 



-aje^ 


( 33 ) 


( 34 ) 


where e = y - y d such that p r + a r . 1 p r ' 1 + ... + a Q e is a Hurwitz polynomial, then the tracking error 
e will converge to zero as time t goes to infinity and the internal state ij will be bounded, || ij || ^ 6, 
where 8 is positive. 

For the multi-input, multi-output system 

x = f(x) ♦ Sj(x) Uj + ... + g m (x)u m 


or 


x = f(x) + G(x)u 


( 35 ) 


y = h(x) 
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For some j-th term j = 1, m 


yj 


- L f h j + 



( 36 ) 


If 


L gi h j(x) = 0 

for all i, we have to differentiate until at least one of the inputs appear in yjj. 
Then 


r ; -l 


y/ = L f J + £ L 8i L f J h j u j 


( 37 ) 


i=i 


In matrix form 


y r = F(x) + B(x) u 


( 38 ) 


y r 



( 39 ) 


where 


F(x) = 


Lf ' h, ... L f rm h r 


( 40 ) 


B(x) = 


r,-l ri -1 

L g , L f h i L g m L f h i 


r - 1 r m -l 

L gl L f h m L g m L f h - 


( 41 ) 


If E‘’(x) exists, control law 


u = -B _1 (x)F(x) + B -1 (x)v 


( 42 ) 
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where 


v 



T 


results in 


y = v 


( 43 ) 


If B _1 (x) does not exist, we must use another method such as a dynamic-equation algebraic method. 


4.3 Continuous Sliding-Mode Control with Linearization 

In classical SMC, chattering is inevitable due to measurement errors, finite switching time, etc. 
In this section, the continuous sliding-mode control with I/O linearization is developed, which may have 
the desirable features of SMC while possibly rejecting the major cause of chatter. For the SISO case the 
output is obtained after differentiating several times until it is related to the input by 
y r = Ljh(x) + LgLj^hfxju. References [21-24] are followed in the approach used here. 

Select a sliding surface 

s»P T Y r _, W 

where 


P T = [Pi - Pr] 
Y r -i -[eft ... e r_1 ] 


s 


P T Y r 


(45) 


and 



T 
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The sliding conditions evolve from 


ss 


= sp T Y r « sp^e ... e r ] 


(46) 


e r = y r 


- y l = L/h + Lf" 1 h ) u - yj 


Using (47) and (46) results in 


ss = sp 


= sp t e 


If 


e ... e^L/h^LgLp 1 h)u-y d r ] 
e + sp 2 e + ... + sP r _! e r_1 + sP r (l/Ii ^gL^" 1 h)u-y d r ] 

= sp r L g L f r_1 h j {u + (sp r L g L f r_1 h) (spj e + sp 2 e + ... sp r _! e r_1 - sp r y d r 

u = -|sp r L g Lf 1 h| (spj e + sp 2 e + ... + sp r _! e r 1 - sp r y d - sp r L g L f h 


(47) 


(48) 


• f t j r- i h ) 2 * o < 49 > 

ss = -[sp r L g L f hj < 0 

For the multi-input, multi-output case, the attractivity condition for SMC towards S = 0, $ = 0 
becomes S T S < 0. 

Let the sliding surface S = [Sj ... s m ] T for some j-th term, j = 1, •••> m 


r j _l 


Sj = E c ji e j‘ 

i=0 


(50) 
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( 51 ) 



r 1 


+ E c j,i-i e / 

i«l 


= Lphj -y d j + Z 

i=l 


r:-l 

Lg L f J h 


T 1 


u i + 


E 


i=l 


In the matrix form 


S = 


r i _1 

E c i,i e i 

i=0 




E C m,i®m 

i=0 


A E = jEj 


- E m ] J 


(52) 


S = e r + E' 


(53) 


where 


e r = 


e l ••• e m 


E = 


r l _1 

E C l,i-l e l 

i=l 


r m -1 


E e m 

i=l 


Using equation (51), for some j (53) results in 


S = L f r h - + E + Bu 


(54) 


where 
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Lf h - y d * E' = 


Lf'^l - »dj * E 


Lf h m y dm + 52 c m,i-l e m 


r l r ml 

y-i - y <\ 


t IT 

E - e; % ... E m 


Ej - E e / > * 1 -- m 


B is the same as equation (41). Therefore, sliding condition 


S T S = E t ILf h - yjj + Ej + E T Bu 


R j ( ET B ) 

B 1 u + -5 — — - 


T„i2 


| E 1 B 


E T lL f r h -y T d +E 


u = E t (l/ h - y d ♦ E'} - k(e t b) 7 


1 E t B | 2 


K is positive definite matrix, then 


s T s = -e t bk(e t b) t < 0 
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4.4 Simulation 


The plant is the simplified second-order longitudinal motion of a high-alpha airplane [1] with 


d = 9.168 c z (a) + q - 1.8366 (8 e + 7) + 7.361904 

(61) 

i| * 5.73 -2.865 

where a is the angle of attack in degrees, q is the pitch rate in degrees per second, and 5 e is the elevator 
control in degrees. 




-1.5 


and 


-0.07378494 a , 

0.09722 a 2 - 2.8653 a + 20.03846 , 
-0.01971 a 2 ♦ 0.74391 a - 7.80753 , 
-0.47333 - 0.01667 a , 


a <£ 14.36 

14.36 <. a <. 15.6 

15.6 £ a ^ 19.6 

19.6 < a < 28 


Let Xj = a, x 2 = q, and the output y = a, then 


X 



= f(x) + g6 e 


(62) 


where 


g T = [-1.8366 -8.595] . 

Uncertainty in f and g were considered as Gaussian noises with variance of unity. 


Here, for convenience and smoothness, we let the desired trajectory 


Yd = “desired = ( A _a o) e " Xt " ( A _ 2 “o) e ^ 


(63) 
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where 


A = 18° 

a c = the equilibrium value with [a D qj = [3.20444632 0] 
S e = -1.8029642 
Xj = 5 , X 2 = -1 
Here, it was assumed that 


e = y - y d 

v = y d + k + I edt + k x e 


(64) 


where k t = -16, k = 10, and I = -164 so that the tracking error will be reduced with a damping ratio 
0.617 and natural frequency 12.96 rad/sec. 

Figures 9 to 14 show the simulated feedback-linearized controlled responses for a change in angle 
of attack a from 0° to 18° in a few seconds according to command (63). These use the control of (64) 
with k = I = 0 for Figs. 9 and 10 and without uncertainty and disturbance for Figs. 9 and 11. Figures 
13 and 14 use SMC with linearization. For SMC with linearization, the sliding surface is given by s = 
Pl e = p t (y-y d ) where p x = 1. Using equation (48) and (62), the control input is u = - spigi- 

The last two figures demonstrate the success of SMC with feedback linearization (nonlinear control) 
relative to the other runs. The constant affine term, derivative, and integral terms are all found to be 
beneficial for both cases. The nonlinear control (NLC) magnitude, however, is the only one within the 
physical constraints when uncertainty and disturbance is considered. 

It must be realized that this is only a first study of a single sliding-mode control used in conjunction 
with feedback linearization. It must be studied with the more complex model for larger variations in 
angle of attack where further dynamic uncertainty is considered along with control rate and acceleration 
constraints. In any event, such controllers can be used for comparison with more robust nonlinear 
controllers and linear controllers. 
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Without 



Figure 9. Feedback-Linearized Control Without Integral Term 





Without I, With disturbance 



Figure 10. Feedback-Linearized Control Without Integral Term With Disturbance 





Line 2 



Figure 11. Feedback-Linearized Control With Integral Term 




With I, With disturbance 



Figure 12. Feedback-Linearized Control With Integral Term and Disturbance 





using SMC 



68P jo O0S/6op 


Figure 13. Sliding-Mode Control (SMC) With Feedback-Linearization (FL) 



using SMC with disturbance 
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Figure 14. SMC-FL Control With Disturbance 



5. OTHER PRELIMINARY STUDIES 
5.1 Bilinear-Based Adaptive Control 

Bilinear systems (i.e., systems linear in state, linear in control but not jointly linear) represent 
perhaps the simplest class of nonlinear systems and, indeed, considers the first nonlinear terms in the 
Taylor series approximation beyond linearization [4], Bilinear systems (BLS) may be further justified 
for aircraft dynamic approximation beyond conventional linearized models since lift and drag aerodynamic 
stability derivatives do depend on control surface variation. For high-performance aircraft using more 
control, such parametric control terms become more dominant, lead to improved controllability and 
performance, and here the BLS becomes a first natural approximation beyond linear systems. Indeed, 
linear analysis fails to predict the desired performance gains and may fail to lead to a stabilizing design 
for sufficiently large rapid maneuvers. While analytical results for nonlinear systems are somewhat rare, 
a considerable base for analysis and design has been developed for BLS [4], 

The MAC design studied in Section 2 uses a generalization of bilinear terms to compute a nonlinear 
model-based adaptive controller. Further analysis will consider the effectiveness of such bilinear -based 
control. But here we consider another approach which could lead to a simple nonlinear feedback 
configuration or could be used as a basis for nonlinear PIF control. 

As a BLS model reference consider the discrete equation 

m 

X k+ 1 = A k x k + B k u k + Y, N ki u ki x k 

i=l 

where, for example, x"^ — [cu q 6 v] is comprised of angle of attack, pitch rate, pitch, and air speed. The 
control u may include, for example, stabilator and thrust vector mechanisms with an error vector taken 
about some desired trim state or trajectory. 

For a first design approach select a quadratic Lyapunov function in ^ (and possibly control rate 
of change). An effective stabilizing controller can then be generated in terms of quadratic and linear 
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feedback terms — the quadratic term being most effective for large errors and the traditional linear 
feedback most effective for small errors [4]. Detailed design and simulation are in progress. 

5.2 Neural-Net-Based Control 

A preliminary study has been started on dynamic aircraft modeling and control by means of 
artificial neural networks. One approach utilizes the time-optimal control study summarized in Section 
3. 

The main idea of the approach is to utilize the approximating capabilities of neural networks to 
extract the closed-loop information from open-loop optimal trajectories. The synthesis of the feedback 
time-optimal controller for a general nonlinear plant is very often untractable from the analytical point 
of view. On the other hand, powerful theory and many efficient numerical algorithms exist for the open- 
loop problem that allow us to calculate optimal trajectories for given initial and final constraints. Those 
trajectories, however, contain information about optimal feedback. For a sufficiently regular optimal 
control problem with integral quality index the principle of optimality assures that an instantaneous value 
of control is a function of the current state, unique up to a set of measure zero. Hence, it may be written 
that 

u ^‘(t) = g(x(t)) 

Each point of an optimal trajectory may be viewed as a point on a hypersurface representing a 
mapping g. If sufficiently many time-optimal trajectories are collected and sampled the mapping g may 
be reconstructed. We propose to use multilayered neural networks to do so. Their capabilities as 
function approximators are widely used in the system and control area to extract complex nonlinear 
functions from discrete measurement data. The theoretical basis for those applications is a result stating 
that given any continuous finite-dimensional function defined over a compact set there exists a feed- 
forward neural network with one hidden layer which approximates the function arbitrarily closely in a 
uniform sense. Since continuous-function compact sets are dense in L 2 , a network approximating optimal 
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control mapping g in a mean-square sense is guaranteed to exist, and one of the training algorithms 
min imizing the mean-square error may be used. The advantage of this approach is that instead of solving 
the nonlinear optimal feedback synthesis problem we deal with an easier problem of calculating open-loop 
trajectories for different initial conditions. 

For the second-order model described above time-optimal control strategies were approximated 
using a simple single-layer neural net. Time-optimal transition was considered to the equilibrium points 
corresponding to angle attack equal 0, 5, 10, 15, and 20. First five separate neural networks were 
trained to approximate the time-optimal strategies for these five target points. To train a network to 
approximate the mapping 

8 ° * = g°»*(a,q) 

a set of open-loop optimizations was performed for different starting points. Then resulting time-optimal 
trajectories were sampled every 0.05 s. For the model in question the time-optimal control is of a bang 
bang type. Training data for the three typical cases is shown for the second-order nonlinear model in 
Fig. 15 with signs corresponding to the control value 0 and " + " signs corresponding to value of -20. 
Figure 15 requires some comment. 

The data was used to train five neural networks with two inputs (a and q), one output (5), and one 
hidden layer of 20 neurons. Sigmoid neurons were used in both hidden and output layer. Each network 
was trained off-line and training lasted 2000 iterations (epochs). In each case, very good approximation 
of function was obtained. 
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Figure 15. Training Data for Three Different Values of Target a with the Set of Equilibria of the System 
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1. INTRODUCTION 

Two nonlinear algorithmic controllers, MAC, are studied here. One uses a 
block-canceling Volterra approximation, and the other MAC consists of solving an 
approximating polynomial time series instate and control. Both methods synthesize 
discrete control sequences and are applied successfully to the control of a simple 
nonlinear longitudinal aircraft model for large variations in angle of attack. 

The Volterra-series approach used here was introduced by Modyaev and Averina 
[1], and a form of inverse generating control according to an assumed structure is 
presented by Harris [2]. This work formed the basis for the methods used here. The 
high angle-of-attack aircraft model derived by Stalford, et al. [3] was the plant 
simulated for the MAC application. In many traditional design studies, a sequence 
of linearized perturbation models are derived for different equilibrium flight 
conditions with linear controllers appropriately derived. Linear adaptive control can 
be derived according to nonlinear gain scheduling of the control law. A highly 
successful version of such control, which includes proportional plus integral plus filter 
(PIF) terms, is presented by Ostroff [4,5], However, such designs usually require a 
large number of set-point design computations, and may have stability problems for 
large fast changes in angle of attack and/or mach number. 

For generation of the nonlinear control, a nonlinear time-series based model 
reference is used. In order to identify such model, experimental data was collected 
for angle of attack (o) and pitch rate (q) subject to random steps of control 


(stabilator, 6). To capture such phenomena as limit cycles in the data the steps were 
rather long (40 s). There were 64 such steps with time discretization of 0.1 s 
resulting in 25,600 points in a state plane for 64 values of control. 

For a least-squares simulated data fit, the following approximation was 

surprisingly accurate: 

o(k+1) = p,„«(k) + P 2a « 2 (k) + P 3a « 3 (k) + 

p 4 „q(k) + p 5 .q(k)a(k) + p 6a q(k)a 2 (k) + p 7 .q(k)a 3 (k) + 
p 8 „u(k) + p 9a u(k)o(k) + p 10 .u(k)« 2 (k) + p 1u u(k)o 3 (k) ♦ p 12a (i) 
q(k+1) = p 1q a(k) + p 2q o 2 (k) + P 3q a 3 (k) + 

P4,q(k) + P 5 q q(k)a(k) + p 6 q q(k)a 2 (k) + p 7 q q(k)a 3 (k) + 
p 8 q qu(k) + P 9 q u(k)a(k) + p 10 q u(k)a 2 (k) + p 11 q u(k)a 3 (k) + p 12q 

Even limit cycles are accurately rendered by this model, as well as the stable zone 
behavior, although large discrepancies occur when the control values are close to the 
stable/unstable zones border. 


2. ADAPTIVE CONTROL APPROACHES 
2.1 Nonlinear Volterra-Based Control 

Here, as in [6], the Volterra series serves as a conceptual starting point for a 
nonlinear time series base control. Continuous time controllers based on Volterra 
series were systematically developed in [7] with formulae for the controller s kernels 
given those of the plant and of the desired feedback system. In particular, the 
problem of so-called exact feedback linearization was solved here. However, those 
formulae may be of limited practical value because of the properties of Volterra 
series under feedback. The problem is that even finite (e.g., second order) Volterra 
series of the open loop results in infinite Volterra series of the closed loop. This 
makes it necessary for the controller to include theoretically an infinite number of 



compensating terms even for a quadratic system. The same problem for the discrete 
time systems was treated in [1] with multidimensional Z transforms to derive the set 
of formulae equivalent to those for so-called exact feedback linearization [8]. 
However, they also provided a very elegant transformation of which results in a 
controller requiring only as many Volterra terms as there are in the assumed plant. 

One attractive feature of this controller is that its structure makes it possible to 
utilize it not only with models represented in the form of Volterra series, but in fact 
with any model with easily divided linear and nonlinear parts of the dynamic 
equations such as (2) above. 

The following algorithm results: 

a) according to the linear part of the plant, calculate the linear control u L (k) 

b) calculate the predicted value of the output at the moment k 

y(k) = L( y(k-1) y(k-M),u(k-lj u(k-M)) 

N{ y(k-1) y(k-M),u(k-1) u(k-M) 

c) solve the "linearizing" control equation for x(k) such that 

N( y(k),y(k- 1 ),...,y(k-M+ 1 ),u L (k) -x(k),u(k-1 ) u(k-M + 1 )) = 

= L(x(k) ,x(k- 1 ) x(k- M + 1 ),9(k),y(k-1 ) y(k- M + 1 )) 

3) calculate the control by 

u(k) = U L (k) - x(k) 

This algorithm becomes a sort of prediction controller which tries to estimate the 
effects of the previous controls knowing the previous values of outputs and then to 
adjust the current value of control so that the nonlinear part of predicted output is 
canceled. 

This discrete time nonlinear a control algorithm is generated according to an off- 
line identification of model (1) with a nonlinear aircraft simulation based on [3], 
Also, a linear controller was designed according to the linear parts of (l)-(3). 

The design was performed to obtain the closed loop model reference behavior 

of the form 


G(z) = 0.05 /(z 2 - 1.6z + 0.65) 

In order not to cancel the zero of the plant, the observer polynomial (z-0.7) was 
introduced. The algorithm for the control value u(k) is as follows. First the estimate 
of the output at moment k is calculated from (1) with k replaced by k-1. 

Then it can be shown that the control becomes 

P 8g U L (k) ~ (P2«° 2 + P3«° 3 + Ps«”P + Pe «^° 2 P 7«^ g3 * Pl2a) (4) 

U(k) + (p 8 « + p 9 .« + Pio«“ 2 + Pn.® 3 ) 

with a(k) and fj(k) designating estimates taken from (1). It is seen that if there are 
no nonlinearities in the model the control reduces to a regular linear controller 

u = u L . 

Simulations were run to test the controller performance especially in the 
unstable range of angle of attack. The system is successfully stabilized and the 
transients are very smooth and without significant overshoots for the nonlinear 
control as demonstrated by Figure la. By different choice of the reference model 
it is possible to obtain much faster, but at the same time much more "nem>us" 
transients. The elevator control is also relatively smooth and within the range 
corresponding to the terminal equilibria. As can be seen from Figure lb, the similar 
linear control is unstable. 


2.2 On-Line Adaptive MAC Algorithm 

i 

Model algorithmic control (MAC), described for example in [2], consists of 
solving the model equation for the value of control necessary to obtain required 
value of output. Usually this desired output trajectory is generated form the setpoint 
by means of a reference model. In case this model is linear, the algorithm in essence 
becomes a linearizing one. 


Here, the controlled output is assumed to be the angle of attack such that the 
reference equation becomes: 

.xjk.1) - S^lk-1) * («(k-1) - a^lk-l)) 151 


with 


*mod(k + D = 

$(k) = [a, a 2 , a 3 , q, qa, qa 2 , qa 3 , U, Uo, Ua 2 , ua 3 , l] T (k) 

a(k) = P« T 4>(k-1) + ( a (k-1) " a mod( k_ ^)) 


q(k) = p q T <t>(k-1) + (q(k-1) - Qmod^- 1 )) 


The controller is assumed to know the values of angle of attack and of pitch rate at 
the moment k-1. Then it estimates their current values o(k) and q(k) taking into 
consideration previous prediction errors and based on them it calculates the control 
required to achieve a ref at the moment k+1. The value of control is found as: 


u(k) = 


a, - p u a - p 2a g 2 - p 3 «a 3 P 4 «q 


p 5 «q“ 


Pe«q« 2 - P/.q® 3 - Pi2« 


Ps« 


+ P 9 .« + Pio .« 2 + Pii «“ 3 


(6) 


where 


« r = a rer (k + 1) - («( k-1) - a^k-1)) 


and a = a(k), q = q(k) as described above. 

The results of the simulations are seen in Figures 2a,b. The-reference trajectory 
was chosen to be l/zM.6z+0.65). The actual output of the plant is seen to follow 
the reference very closely, even though the region of operation was that of the most 
severe nonlinearities. The control action is also remarkably smooth. 

The discrete time nonlinear state space model (1) describes the behavior of the 
complex nonlinear plant quite accurately in the entire region of operation. In 
practice, however, such a global model is rather difficult to fit, and consequently one 


should look for local approximations, depending on the current operating conditions. 
In such a situation, on-line adaptive control seems to offer an ideal solution. 

The algorithm discussed in the previous section can be made adaptive, or self- 
tuning, by incorporating on-line identification of the parameters. A recursive least 
squares (RLS) algorithm was implemented in the following form taken from [8]. 

Dfl0 e(k) (7) 

X( k-1) + 4>(k-1) T Q(k-2)<j>(k-2)4»(k-1) 


Q(k-1) 


1 

Mk-1) 


Q(k-2) 


Q(k-2) <}>(k-1)<t>(k-1) T Q(k-2) ) 
A(k-1) + 4»(k- 1 ) T Q(k-2) 4>(k-1) J 


( 8 ) 


©(k-1) = y(k) - p T «l>(k-1) 


where y may denote a or q and p may stand for p. or p q , respectively. The forgetting 
factor X was introduced to enable the algorithm to change the estimates of 
parameters with the change of operating conditions. To avoid the unlimited growth 
of covariance matrix Q at the steady state when the input is not persistently exciting 
the variable forgetting factor policy was implemented: 


m = i 


e(k) 2 

e(k) 2 


( 10 ) 


where e(k) is the current prediction error, e(k) is the average prediction error form 
last 10 samples and e is equal to 0.01. As an additional precaution the trace of the 
covariance matrix Q was monitored and Q was reset to diagonal matrix whenever the 
threshold value was exceeded. Starting values of parameters were taken to be as in 

( 1 ). 

Figure 2 displays the simulation results for a reference model specified as 
l/(zM.8z+0.82). Remarkably exact following of the reference trajectory may be 
observed, although, surprisingly enough, the performance is slightly worse than in the 
nonadaptive case. Most probably this is due to the fact that prediction error now 
changes much more quickly because of the ongoing identification process. Thus, 
approximating the term (y(k+ l)- y<mxl (k+ 1)) by (y(k-l)-y_(k-l)) may worsen the 


behavior of the system as two values of no longer correspond to the same 
parameter vector. Since the on-line identification process assures (at least in 
principle) that the prediction error should asymptotically converge to zero it is 
possible that the correction terms in a(k), q(k), and in control equation (5) ought to 
be omitted. 

The performance of the adaptive nonlinear MAC controller was compared to 
the linear one, which uses the same control strategy but with a strictly linear model 
being identified and used for the calculation of the control action. Clear difference 
between the performance of linear and nonlinear controller can be seen from Figure 
3, particularly in control action at the setpoint a = 15°. The linear identifier has 
obvious difficulties with fitting the parameters of a linear model to the behavior of 
the plant which is highly nonlinear in this region. As a result, the control starts 
oscillating for a while. Also, it was seen that the nonlinear algorithm results in 
control plots that are more smooth, although they still contain one-pulse spikes. To 
eliminate these spikes weighting of the increments of control can be introduced into 
the algorithm with little performance deformation. 


4. CONCLUSIONS 


The nonlinear control applications to high angle-of-attack aircraft, as reported 
here, is of a preliminary nature. However, the analysis does suggest that nonlinear 
adaptive control can be quite effective to stabilize large rapid maneuvers in angle of 
attack. Of the comparisons made, the on-linear, nonlinear-time-series and adaptation 
performed the best and was quite superior to a similar linear MAC. 
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Summary 

The stability analysis of an airplane using its nonlinear model is 
presented. The analysis is based on the robust stability analysis 
approach for linear systems. Then, based on analysis, a small 
static feedback gain is designed such that the robustness of the 
closed-loop nonlinear system stability is significantly improved. 
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1 . INTRODUCTION 


The stability is one of the most important issues in the control 
system design. Recently there has been observed a great interest 
in the methodology of robust stability analysis and design of 
robust control systems for linear dynamic systems [6] . The 
objective of this paper is to investigate the applicability of 
this approach for nonlinear dynamic system such as an aircraft 
flight in high angle of attack/sideslip flight. The unstable 
control system can result in the plane crash. 

There is considered stability of nonlinear, simplified however, 
model of the airplane. The organization of the paper is as 
follows. In section 2, the model of the plane is presented. 
Stability of the aircraft is considered in section 3. Finally, 
concluding remarks are given . 


2 . THE AIRCRAFT MODEL . 


Model of 
usually, 
design, 
given in 


an airplane is highly nonlinear, [4,5]. There are 
however, used simplified models for control system 
e.g. [1,3,9]. In this paper we consider very simple model 


[ 8 ] : 


x=A (x) x+Bu+D 


( 1 ) 


where x= 
degrees , 
elevator 


is a state vector, 

qj 

'q is the pitch ratio m 
control in degrees. 


a is the angle of attack 
degrees per second and u is 


in 

the 


9.168c (a) 1 

-5.73 0 

and c (a) is a nonlinear function. This function can, however, be 
z 

approximated as follows: 


_ [-1.83361 [-5. 473296" 

B= [-8 . 5 9 5 0 J ' U "L 2 . 865000_ 


2 


c z (a) 


f-0 .072815870 

0.088470922-2.3774/a 
‘ 0.033099050-1.4068/a 

-0 . 016633734-0 .4743 /a 


for o; 
for 14.74° 
for 17.40° 
for 18.87° 


< a s 14.74° 

< a s 17.40° 

< a s 18.87 

< a s 28.00° 


( 2 ) 


It is easy to find that 


-0.072818087 < c z (a) * -0.048161261 
-0.047751524 < c z (a) * -0.041453149 
-0.041768869 < c z (a) * -0.033573019 


for 

14.74° 

< 

a 

< 

17.40 

for 

17.40° 

< 

a 

< 

18.87 

for 

t— * 
00 

♦ 

00 

o 

< 

a 

< 

28.00 


This model approximates model taken from measured wind tunnel 
values of the T-2C airplane [7]. It is known that numerical values 

of c (a) and b~ are uncertain, 
z z 


Our purpose is to consider stability of system (1) in the range of 
angle of attack 0°sa£28°, and to find a static feedback which can, 
eventually, improve the stability of the plane m this range. 


3. STABILITY ANALYSIS. 

Consider linear time-invariant system 


where xeR n is a state vector. Then, assuming that the system is 
asymptotically stable one can define the following notions, [2] . 

Definition 1 . 

A connected set Dj in the system parameters-space (parameters o 
matrix A) is a robust time invariant stable (RTIS) set £or system 
(4) iff AeD-j. and every t ime- invariant system 


is asymptotically stable for jtfeDj . 

D 
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Definition 2 . 


A connected set By in the system parameters -space is a robust time 
varying stable (RTVS) set for system (4) iff Aefly and every time- 
varying system (5) is asymptotically stable for 


Then, consider four linear models instead of (2), respectively: 


-0.072815870 

for 

-0.048161261 

for 

-0 . 041453149 

for 

-0.016633734 

for 


0° 

3 

a 

< 

14.74 

14.74° 

< 

a 

< 

17.40 

17.40° 

< 

a 

< 

18.87 

18.87° 

< 

a 

< 

28.00 


It is easy to find that all models are asymptotically stable. We 
are, however, interested in the set of (k 1# k 2 ) such that all the 
linear closed loop systems will be stable with the following 

feedback 

u=Kx, K=[k^ k 2 ] 

An appropriate region ^ can be easily calculated based on 
algorithm 2 proposed in [2]. This is, however, only the second 
order system and one can simply obtain analytical formulas for the 
RTIS region in this case. The characteristic polynomial for the 
2nd order system has the following form 

s^+as+b=0 


It is known that all roots of this polynomial are in the left half 
plane, i.e. a system is asymptotically stable stable, iff a>0 and 
b>0 . Based on this, the RTIS region n z for 'stable' feedback gains 
was calculated. The region is presented on fig. 1, a dashed line 
represents RTIS region for model PI, 0°sasl4.74°, a dotted model 
P2 for 14 . 74 °<a^l7 .4°, a dash-dotted model P3 for 17.40 <asl8.87 
and a continuous line model P4 for 18.87 <a^28 . It is easy to see 
that the system without feedback, i.e. k 1 =k 2 =0, sign + on the 
plane (k-,^), is very close to the stability region boundary. One 
can improve stability assuming appropriate K from ^(kyk^. 
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Next, RTVS sets £1^ were calculated for these models, according to 
the algorithm given in [2], for uncertain parameters a n and a 21 
in A. They are presented on fig. 2. All four models are inside the 
RTVS region calculated for the model P4. Moreover, since all time 
varying (nonlinear) a n = 9 . 168c z (a) is smaller than nominal values 
used in linear models it means that the whole nonlinear system (1) 
is asymptotically stable for 0°sas28°. However, there is a very 
small upper bound for a^ in this model, namely 

+A a H < 0.0227 

This can cause that with small system uncertainty the system can 
be unstable. The vertexes of RTVS quadrilateral Qy on the plane 
(Aa 11 ,Aa 21 ) are as follows: 

= { (-176.5,0), (0.0226,0), (0,-0.2274) (0,0.2944) } 

In order to improve system stability feedback gain matrix K was 
chosen from Oj. Intuitively, it seems that a good gain is a small 
one - a high gain can result in a lack of system controllability 
because of saturation of the control input, and such that the 
stability margin with respect to K will be rather large. 

Thus, the good choice seems to be K^=[0 0.2] . For this gain one 
obtains significant improvement of RTVS set. This set is shown on 
fig. 3. In this case also all models are inside the £1^ calculated 
for the model P4 . However, an upper bound for a^ is more than 8 
times greater: 

+Aa,. < 0.1835 

Also range for uncertain parameter a 21 is almost 6 times larger. 
Indeed, the vertexes of RTVS quadrilateral in this case on the 
plane (Aa^Aa^) are as follows 

V = { (-25049,0), (0.1835,0), (0,-3.547), (0,3.211) } 
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Then, it was considered feedback gain 0.2 0]. This gam, 

however, seems to be worse situated in the RTIS set ^ than K 1 
considering the stability region with respect to K. Nevertheless, 
also in this case one obtains improvement of robust stability for 
the closed loop system. The appropriate RTVS set ^ is shown on 
fig. 4. In this case an upper bound for a^ is as follows 

+ ^ a ll < 0-0613 

Similarly, range for perturbation in a 21 is larger than for K=0. 
The vertexes of RTVS set are as follows 

V = { (-65.15,0), (0.0613,0), (0,-0.5770), (0,1.3641) } 

V2 

It should be noted that all RTVS sets were calculated under 
assumption Q=I in algorithm 3 [2] . 

From the above analysis follows that relatively small static 
linear feedback gain K=[0 0.2] significantly improves stability 

of the system. It should be emphsized that every nonlmear/time- 
varying system (1) with a n and a 21 from the obtained RTVS set 0^ 
will be asymptotically stable. This way we have designed a robust- 
stable nonlinear closed-loop system. 


4 . CONCLUDING REMARKS . 

A robust-stable nonlinear controlsystem has been designed. It is 
shown that small linear static feedback gain can significantly 
improve stability of the airplane. The feedback gain seems to be 
so small that it should not constrain control signal during plane 
maneuvering. This should also results in better a controllability 

of the plane. 

The stability analysis and feedback gain synthesis were done using 
methods designed for linear systems [2]. This, approach can be 
also used for more complicated nonlinear systems. For instance, 
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assuming as a base model for the airplane, the linear 9th order 
model given in [1,9]. This model is unstable, but, as it was shown 
in [2], one can deal also with unstable models using the same 
approach . 

Presented results also show the power of the approach proposed in 
[ 2 ] . 
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